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Adsorption of para-hydrogen on the surface of graphite pre-plated with a sin- 
gle atomic layer of krypton, is studied theoretically by means of path integral 
ground state Monte Carlo simulations. We compute energetics and density 
profiles of para- hydrogen, and determine the structure of the adsorbed film 
for various coverages. Results show that there are two thermodynamically 
stable monolayer phases of P-H2, both solid; one is commensurate with the 
krypton layer, and has coverage 6 « 0.0636 A" 2 ; the other is incommensu- 
rate, has coverage 9 ~ 0.0716 A~~ 2 , and is compressible up to 6 ~ 0.0769 
A~ 2 . No evidence is observed of a liquid phase at T=0 for intermediate cov- 
erages. These results are qualitatively similar to what is seen for p-B.2 on 
bare graphite. Quantum exchanges of hydrogen molecules are suppressed in 
this system. 

PACS numbers: 02.70.Ss,67.40.Db, 67.70. +n,68.43.-h. 



1. INTRODUCTION 

Low temperature adsorption of highly quantal fluids, such as helium or 
para-hydrogen (p-H^), on a variety of substrates, has been the focus of much 
experimental and theoretical work. A major motivation of these studies is 
the exploration of the fascinating properties that such adsorbed quantum 
films display, often considerably different than those of the bulk materials. 

In particular, a fluid of P-H2 molecules is an interesting physical system 
for a number of reasons. Because a P-H2 molecule has half the mass of 
a helium atom, zero-point motion can be expected to be quite significant; 
each molecule is a spin-zero boson, and therefore it is conceivable that, at 
low enough temperature, a P-H2 fluid might display physical behavior similar 
to that of fluid helium, including superfluidity.^ 
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Unlike helium, though, bulk P-H2 solidifies at low temperature (T c ~ 
14 K); this prevents the observation of phenomena such as Bose Einstein 
condensation (BEC) and, possibly, superfluidity (SF), which are speculated 
to occur in the liquid phase below T ~ 6 K. Solidification is due the depth 
of the attractive well of the potential between two hydrogen molecules, sig- 
nificantly greater than that between two helium atoms. Several, attempts 
have been 

made pEEl 

to supercool bulk liquid P-H2, but the search for SF 
(in the bulk) has so far not met with success. 

Potential avenues to explore, toward stabilizing a liquid phase of P-H2 
to low enough temperatures that a SF transition could be observed, include 
the reduction of dimensionality. An extensive theoretical study of the phase 
diagram of P-H2 in two dimensions (2D) has been recently carried out, based 
on Path Integral Monte Carlo (PIMC) simulations.^ The main result is that 
the equilibrium phase of the system at low T is a triangular crystal, with a 
melting temperature T m ~ 6.8 K. This value is approximately half that of 
bulk p-H-2, but still significantly higher than the temperature at which the 
system, if it remained a liquid, would undergo Bose Condensation and turn 
superfluid, estimated at ~ 2 K in 2D. Another important result of the same 
study, is that no metastable liquid phase exists; for, the system remains a 
liquid down to the spinodal density, below which it breaks down into solid 
clusters. These results, in part, cast doubts on the prospects of observing a 
liquid (SF) phase of P-H2 at low T, in 2D. 

The closest experimental realization of a two-dimensional system, is a 
film of P-H2 molecules adsorbed on a substrate. However, quantum zero- 
point motion of adsorbed particles in the direction perpendicular to the 
substrate can be quite significant, as calculations for adsorbed helium films 
on alkali metal substrates showP One might speculate that such zero-point 
motion may result in an effective screening of the P-H2 intermolecular in- 
teractions, possibly leading to a reduction of T m , with respect to the purely 
two-dimensional case. Indeed, some PIMC studies^ of P-H2 films adsorbed 
on a lithium substrate yielded a melting temperature of approximately 6.5 
K, i.e., slightly lower than in 2D. The interesting physical question is whether 
T m may be reduced even further, by an appropriate choice of substrate, to 
the point where collective quantum many-body phenomena could become 
observable, in the liquid phase. 

Experimentally, adsorbed films of P-H2 on surfaces are readily acces- 
sible, and indeed have been studied extensively on a variety of substrates. 
For example, various have been used to study the phase 

diagram and structure of monolayer p-U.2 films adsorbed on graphite. In a 
recent neutron scattering investigation of deuterium (D2) films adsorbed on 
a krypton pre-plated graphite substratepl evidence of a stable liquid phase 
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of D2 down to T ~1.5 K was reported. Such a result is truly remarkable, 
particularly considering that it pertains to the heavier isotope of hydrogen, 
which should display an even stronger tendency to crystallize than P-H2. 
Motivated by this finding, we have undertaken a theoretical study of the low 
temperature phase diagram of P-H2 films adsorbed on krypton pre-plated 
graphite. 

On general grounds, one would expect p-B.2 to remain liquid as well, 
on account of its lighter mass, even though the inter molecular interaction is 
slightly more attractive for p-H^U^^We explore a fairly simple model of our 
system of interest in its ground state (T=0 limit) to serve as a baseline for fu- 
ture, more refined calculations. Energetic and structural properties of a layer 
of P-H2 molecules adsorbed on Kr /Graphite are investigated theoretically, 
by means of Path Integral ground state (PIGS) Monte Carlo simulations. 

The main results of this study are the following: 

1. A stable commensurate solid P-H2 monolayer exists, with coverage (i.e., 
2D density) # o =0.0636 A -2 . This solid is commensurate with the kryp- 
ton plating. 

2. An incommensurate solid P-H2 monolayer exists, with coverage 9\ =0.0716 
A~ 2 , compressible up to # 2 =0.0769 A" 2 . 

3. No evidence is observed of a thermodynamically stable liquid phase in 
the ground state of this system. 

4. Evidence of quantum exchange between p-Yi.2 molecules is not found 
in the T=0 limit; this in turn indicates absence of superfluid behavior 
in this system. 

The remainder of this manuscript is organized as follows: Sec. |2~1 offers 
a description of the model used for our system of interest, including a discus- 
sion of the potentials and justifications for the main underlying assumptions. 
Sec. l3~l involves a brief discussion of the computational technique and specific 
details of its implementation, in addition to details of calibration and opti- 
mization. The results are presented in Sec. El finally, Sec. ED is a summary 
of the findings and our concluding remarks. 

2. MODEL 

We consider a system of N P-H2 molecules, sitting above a substrate 
consisting of a single atomic layer of krypton, below which is a graphite 
substrate, which we assume to be smooth. The (L) krypton atoms are as- 
sumed fixed in space at positions R^, owing to their relatively large mass. 
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They sit at a height of 3.46 A over the graphite; this distance corresponds 
to the minimum of the most accurate Kr-graphite interaction potential cur- 
rently available! 16 * 17 * All of the Kr atoms and P-H2 molecules are regarded 
as point particles. The model quantum many-body Hamiltonian is therefore 
as follows: 

N N L N 

h = -2^ E v ? + E v(r*) + EE^- R ^l)+E u( Zi ) (1) 

i=l i<j i=l k=l i=l 

Here, m is the mass of a P-H2 molecule, {r,} (with j=l,2,...,N) are the 
positions of the P-H2 molecules, = |r» — |; Zi is the height of the ith 
P-H2 molecule above the graphite surface. V is the potential describing the 
interaction between any two P-H2 molecules, and U represents the interaction 
of a P-H2 molecule with a Kr atom. 

As mentioned above, the underlying graphite substrate is regarded as 
smooth. The justification for this assumption is that, due to the presence 
of the Kr spacer layer, the P-H2 molecules remain at a distance of at least 
~ 7 A from the graphite substrate. Thus, the effect of the corrugation of 
the substrate should be negligible. Therefore, we use a simple "3-9" poten- 
tial to describe the interaction of p-H.2 molecules with the smooth graphite 
substrate.™ U has the following form: 

= - 7 <*> 

where C=7913.24 A 3 K and £>=259.39 K are parameters derived from the 
original p-H-2-C Lennard-Jones parameters^ (<j = 3. 18 A, e = 32.05 K) and 
the density of carbon atoms in graphite (p = 0.114 A~ 3 ). 

All pair potentials are assumed to depend only on relative distances. 
The interaction V is described by the Silvera-Goldman potentialpD which 
provides an accurate description of energetic and structural properties of 
condensed P-H2 at ordinary conditions of temperature and pressure! 22 * 23 * 

The interaction of a p-H-2 molecule and a Kr atom is modeled using a 
standard 6-12 Lennard-Jones (LJ) potential; to our knowledge, there are no 
published values for the parameters of this potential for this specific inter- 
action. Therefore, we make use of the Lorentz-Berthelot mixing rule, 24 25 
yielding e = 75.6 K and a = 3.3 A. 

The model (^Q) clearly contains important physical simplifications, such 
as the neglect of zero-point motion of Kr atoms, as well as the restrictions to 
additive pairwise interactions (to the exclusion of, for example, three-body 
terms), all taken to be central, and the use of the highly simplified LJ and 
"3-9" potentials. Nonetheless, it seems a reasonable starting point, and even 
quantitatively we expect it to capture the bulk of the physical picture. 
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3. COMPUTATIONAL METHOD 

Accurate ground state expectation values for quantum many-body sys- 
tems described by a Hamiltonian such as Q can be computed by means 
of Quantum MOnte Carlo (QMC) simulations. In this work, the method 
utilized is Path Integral ground state (PIGS), which is an extension to zero 
temperature of the standard, Path Integral Monte Carlo method!^ PIGS is 
a projection technique, which filters the exact ground state wave function 
out of an initial trial state. It is therefore closely related to other ground 
state projection methods, such as Diffusion Monte Carlo (DMC), but has a 
few distinct advantages (for a discussion, see, for instance, Ref. I27)) . Because 
this method is described in a number of publications, it will not be reviewed 
here. Some of the technical details of the calculation performed in this work 
(mainly, the short imaginary time propagator) are the same as in Ref. 1281 

The trial wave function utilized is of the Jastrow type: 

,N \ / N L \ / N \ 

^Hri,r 2 ,...r^)= (ne-^Mx(nn e ^ (|r ^ Rfcl M x ( ( 3 ) 

\<j ' \=\ k=l ' \=\ ' 

with pseudo-potentials w (p-H2-graphite), u (p-H2-Kr), and v (P-H2-P-H2) 
chosen as follows: 

w (r) = fL , u (r) = 1 and v(r) = ^ (4) 

The values of the parameters a = 30 A 3 , 7 = 250 A 5 and /x = 750 A 5 
were obtained empirically, by minimizing the energy expectation value com- 
puted in separate variational calculations. Using the trial wavefunction as 
defined above, we observe convergence of the ground state energy estimates 
with a projection time 0.250 K _1 , using a time step r = 7.8125 x 10 -4 K _1 . 

PIGS calculations for a range of P-H2 coverages were carried out, start- 
ing from an initial configuration of para-hydrogen molecules sitting atop 
the Kr layer. The simulation cell consists of a 6x8 triangular lattice of Kr 
atoms with 4.26 A nearest neighbor spacing. Periodic boundary conditions 
are used in the three directions, but the simulation cell is chosen sufficiently 
big in the z direction that they have no effect vertically. Because of the 
strongly attractive character of the composite substrate, for a small enough 
number of hydrogen molecules, the system remains vertically confined (i.e., 
P-H2 molecules do not evaporate). 

The systematic errors of our calculation are attributable to finite pro- 
jection time and the finite time step r. Based on comparisons of results 
obtained from simulations with different values of projection time and time 
step, we estimate our combined systematic error on the total energy per 
P-H2 molecule to be of the order of 0.7 K or less (below 0.5%). 
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Fig. 1. Energy per P-H2 molecule e(N) (in K) computed by PIGS, as a 
function of the coverage 6 (in A -2 ). 



4. RESULTS 

Physical quantities of interest include the ground state energy per P-H2 
molecule, e(N), and the vertical P-H2 density profile, p(z), above the com- 
posite substrate. These quantities can be computed in an unbiased fashion, 
i.e., the variational bias arising from the initial state (i.e., trial wave function 
P)l) can be essentially completely removed.^ 

Results for e(N) are shown in Fig. ^ The main features are two 
energetic minima, one at o = 0.0636 A -2 , a coverage corresponding to 
commensuration (iV=48 P-H2 molecules). The second, at 9\ = 0.0716 A~ 2 , 
corresponds to an incommensurate solid monolayer, which remains stable, 
based on an analysis of the associated chemical potential (see below), up to 
2 = 0.0769 A~ 2 (iV=58 p-B 2 molecules). 

A direct comparison can be made between our e{N) curve, and that ob- 
tained by Nho and Manousakis (Ref. 29 ), who studied P-H2 monolayer films 
on bare graphite (at low temperature). The shape of the energy curve in both 
cases is very similar; there is an energetic minimum at commensuration, cor- 
responding to precisely the same coverage, followed by a negative-curvature 
increase. Where we find a second energetic minimum, they find a change in 
curvature, in both cases preceding a positive curvature (thermodynamically 
stable) portion, although, on Kr pre-plated graphite, such a region seems 
more extended. 

We calculate the chemical potential, shown in Fig. [21 by first fitting 
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potential, /U (in K), computed as explained in the text. Dotted vertical line 
shows the chemical potential of bulk solid P-H2, at the T=0 equilibrium 
density (from Ref. I23]) . 



the curve for e(N) and then minimizing the grand canonical energy 4>(N) = 
N(e(N) — /i) with respect to N, for different values of /u. As shown by the 
data in Fig. [21 these are the only stable coverages, at least up to 6 « 0.0848 
A~ 2 (the highest coverage explored in this work). The presumption is that 
the next thermodynamically stable configuration would be at second layer 
completion, 6 w 0.127 A" 2 . 

Fig. |31 shows the (three-dimensional) p-U.2 density profile p(z) for a 
coverage O (commensurate solid layer), as well as the P-H2 density profile 
for the highest stable coverage, 62 (incommensurate solid layer). Also shown 
in Fig. |31 is the density profile for p-^2 on lithiurrP at low temperature (2 
K). It is evident that for the substrate studied here, P-H2 is more localized in 
the direction perpendicular to the substrate (the density profiles are peaked 
much more strongly) in comparison, despite being at T=0. Thus, the physics 
of p-13-2 is considerably more 2D on this substrate, than on the weak Li 
substrate. This is consistent with the stronger substrate attraction, and 
leads us to predict a melting temperature for the adsorbed film of the order 
of ~ 7 K, i.e., comparable to that of purely 2D p-H-2- Obviously, however, 
only finite temperature calculation can address this point quantitatively. We 
also see that the profile for 62 is extended slightly further to the right than 
that for O ; in accommodating more P-H2 in the first layer, the molecules are 
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Fig. 3. Density profile p(z) (in A~ 3 ) of adsorbed P-H2 at a coverage 9 (solid 
line with lower maximum value), and at 82 (solid line with higher maximum 
value). The square represents the position of the graphite substrate and the 
circle that of the Kr layer. Dashed line is the density profile for P-H2 on 
a lithium substrate at the coverage 0=0.070 k~ 2 (from Ref. El it has been 
shifted along z for direct comparison). Distances are expressed in A. 



"squeezed" to occupy regions higher above the substrate. 

Snapshots of typical configurations of the P-H2 solid film are displayed 
in Figs. ^ and [5J The first, at O , shows the commensurate structure of 
p-H-2 corresponding to this system's first energetic minimum. The second is 
a snapshot at the highest stable coverage, 62', P-H2 form an incommensurate 
solid, not found to be rotated relative to the Kr lattice (although the size of 
our simulated system clearly limits our capability of resolving such an issue) . 
56 of the 58 molecules in the simulations appear to try to form a 7x8 regular 
triangular lattice, with the two additional molecules packing in, giving rise 
to the clear dislocations. Simulations with higher coverages (thermodynam- 
ically unstable) were consistently found to put additional particles in the 
second layer. 

In contrast, for P-H2 on bare graphite,^ the somewhat stronger ad- 
sorption potential allows for an even denser packing; compression of the 
incommensurate phase is reported well beyond 9 = 0.0849 A" 2 , with the p- 
H2 lattice rotated 5° relative to the graphite lattice. We should also mention 
that, due to the system size and geometries employed in this work, we could 
not observe such crystal phases as the uniaxially compressed and the stripe 
one, which have been experimentally observed and theoretically studied for 
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Fig. 4. Snapshot of a typical configuration of p-rFj molecules adsorbed to 
the Kr/graphite substrate at a coverage of 6> o =0.0636 A~ 2 . Large circles are 
Kr atoms. The positions of all molecules at each one of 320 imaginary time 
slices are shown as discrete paths. Distances are expressed in A. 
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Fig. 5. Same as Fig. 0]but at a coverage of 6*2=0.0769 A" 2 
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p-H-2 on graphite)^ we presume that these phases should exist in this sys- 
tem as well, and can certainly be studied with the computational method 
used here. However, based on the results of Ref. ,221 we do not think that 
their inclusion would significantly alter our main conclusion, concerning the 
existence of a liquid phase, which was our main interest in this work. 

It should also be noted that the computational method adopted here 
does not allow one to make a direct estimation of the P-H2 molecules ex- 
change frequency, unlike its finite temperature counterpart (Path Integral 
Monte Carlo). Nevertheless, visual inspection of many-particle configura- 
tions generated in the Monte Carlo simulation shows little or no overlap of 
paths associated to different molecules, which is substantial evidence that 
many-particle permutations are absent in this system. This is consistent 
with the high degree of localization that molecules experience. 



5. CONCLUSIONS 

Using a numerically exact ground state Quantum Monte Carlo method, 
we studied P-H2 adsorption to krypton pre-plated graphite. We performed 
calculations based on a simple model, in which graphite corrugation is ig- 
nored, the Kr atoms in the spacer layer are assumed static and point-like, 
and p-H2-substrate interactions are given by Lennard-Jones type potentials. 

We find that there are two stable phases of P-H2, both solid; one is a 
monolayer commensurate with the Kr layer, while the other is an incommen- 
surate monolayer, compressible within a small range of coverages. Quantum 
exchanges of hydrogen molecules are suppressed in this system. This is sim- 
ilar to what is seen for hydrogen on bare graphite. 29 Altogether, this study 
has provided no evidence of a thermodynamically stable liquid phase of P-H2 
at T =0 on the substrate considered here; it is unlikely, based on compar- 
isons with P-H2 on lithium from Ref. |S| and 2D P-H2 from Ref. H3 that a 
liquid would appear, in our model, at temperatures as low as 1.5 K. 

There are obviously several sources of uncertainty in this calculation 
which need to be discussed. The potentials used to describe the interactions 
between the p-H-2 and the substrate are very rough; this does not seem 
too important an issue, as far as the interaction of P-H2 molecules with 
graphite is concerned, given the relatively large average distance at which 
molecules sit. On the other hand, a more realistic interaction potential 
between P-H2 and krypton may conceivably alter the basic energetics shown 
here. Despite these issues, and other simplifications, it does not seem likely 
(to us) that the basic structural information will change dramatically. It is 
also unlikely that the basic physics will change much at temperatures as low 
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as ~ 1.5 K (namely those reached by the experiment of Ref. I13|) . though 
finite temperature calculations are planned to address this concern, as are 
simulations of larger systems. 

Thus, our preliminary conclusion is that our calculation appears to yield 
results in disagreement with the experimental findings of Ref. ..dj reporting 
a low temperature liquid phase for the heavier isotope of P-H2 (D2), which 
should display an even stronger tendency to crystallize than p-B.2 , on account 
of its heavier mass. 
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